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Abstract 

We  derive  a  control  theoretic  model  of  sector-based  air  traffic  flow  using  hybrid 
automata  theory.  This  model  is  Lagrangian,  meaning  that  it  models  the  properties  of 
the  system  along  its  trajectories.  A  subset  of  this  model  is  used  to  generate  analytic 
predictions  of  air  traffic  congestion:  we  define  and  derive  a  dynamic  sector  capacity 
which  we  use  to  predict  the  time  it  takes  to  overload  a  given  portion  of  airspace.  This 
result  links  our  approach  with  Eulerian  models,  which  account  for  temporal  variations 
of  parameters  in  a  fixed  volume.  We  design  and  validate  an  air  traffic  flow  simulator, 
to  assess  the  accuracy  of  our  predictions.  The  simulator  is  then  used  to  show  that  flow 
scheduling  and  conflict  resolution  may  be  decorrelated  under  assumptions  on  aircraft 
density. 
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Symbol 

Definition 

0% 

initial  arc  length  distance  of  aircraft  i  from  the  SFO  airport  along  arrival  route 

a 

bi 

vector  [a?,  •  •  •  ,  a^]  of  initial  arc  length  distances  for  the  N  aircraft 
variable  used  to  compute  the  mode  switching  time  of  aircraft  i 

b 

vector  ,  6jv]  for  the  N  aircraft 

(p  . 

mm 

minimal  distance  from  aircraft  i  to  any  other  aircraft  in  the  sector 

distLOs 

distance  at  which  a  LOS  occurs 

A  L 

requested  outflow  separation  for  merging  traffic  in  the  region  of  interest 

AL;n 

imposed  inflow  separation  for  cross  traffic  in  the  region  of  interest 

A  Lout 

requested  outflow  separation  for  cross  traffic  in  the  region  of  interest 

A  Tout 

outflow  rate  of  the  region  of  interest  (one  aircraft  every  ATout  time  units) 

ATact 

time  period  for  one  iteration  of  the  simulator 

^^LOS 

time  until  the  next  predicted  occurrence  of  a  loss  of  separation  for  aircraft  i 

/(•) 

penalty  function  for  aircraft  separation  (associated  to  Lmin  for  all  i) 

J 

cost  function  encoding  ATC  controller  action  and  sector  state 

given  action} 

cost  associated  to  a  (given  action}  of  the  controller  (VFS,  shortcut,  etc...) 

M 

Mach  number 

N 

1  v  max 

total  number  of  aircraft  in  the  sector  of  interest  at  a  given  time 

^maneuver 

number  of  maneuvers  the  simulator  can  assign  to  an  aircraft  at  any  given  time 

-^choice 

number  of  aircraft  selected  for  analysis  by  the  simulator 

-Amoved 

total  number  of  aircraft  moved  at  a  simulator  iteration 

nLOS 

number  of  LOS  for  aircraft  i  which  would  happen  with  a  given  set  of  maneuvers 

^limit 

dynamic  capacity  of  the  sector  of  interest 

rotation  matrix  of  angle  ijj  for  heading  changes 

rpi 

^  breach 

boundary  condition  breach  time  of  aircraft  i 

./.switch 

Li 

time  at  which  aircraft  i  undergoes  a  mode  switch  by  ATC 

t  block 

time  at  which  a  metering  condition  is  imposed 

rr 

-*■  limit 

saturation  time  of  the  sector  of  interest 

TOAlred 

predicted  time  of  arrival  of  aircraft  i  (at  TRACON) 

TOAl eal 

actual  time  of  arrival  of  aircraft  i  (at  TRACON) 

^current  heading 

velocity  vector  of  a  given  aircraft  at  its  current  heading 

^min 

minimum  aircraft  speed 

^nom 

nominal  aircraft  speed 

^max 

maximum  aircraft  speed 

^given  action 

weight  (penalty)  associated  to  a  given  action  of  the  controller 

Xi 

position  of  aircraft  i 

Xex 

end  of  the  controlled  area 

X 

distance  to  the  destination  airport 

^  switch 
xi 

location  at  which  aircraft  i  undergoes  a  mode  switch  by  ATC 
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Acronym 

Definition 

ARTCC 

Air  Route  Traffic  Control  Center 

ATC 

Air  Traffic  Control 

ATM 

Air  Traffic  Management 

ATCSCC 

Air  Traffic  Control  System  Command  Center 

BC 

Boundary  Conditions 

CONUS 

Continental  US 

ETMS 

Enhanced  Traffic  Management  System 

FACET 

Future  ATM  Concepts  Evaluation  Tool 

LWR 

Lighthill-Whitham- Richard 

LOS 

Loss  Of  Separation 

NAS 

National  Airspace  System 

nm 

Nautical  Mile 

SFO 

San  Francisco  Airport 

TFM 

Traffic  Flow  Management 

TMU 

Traffic  Management  Unit 

TRACON 

Terminal  Radar  Approach  CONtrol 

VFS 

Vector  For  Spacing 

1  Introduction 

The  National  Airspace  System  (NAS)  is  a  large  scale,  layered,  nonlinear  dynamic  system: 
its  control  authority  is  currently  organized  hierarchically  with  a  single  Air  Traffic  Control 
System  Command  Center  (ATCSCC),  in  Herndon  VA,  supervising  the  overall  traffic  flow. 
This  is  supported  by  22  (20  in  the  continental  US  or  CONUS)  Air  Route  Traffic  Control 
Centers  (ARTCCs,  or  simply,  Centers)  organized  by  geographical  region  up  to  60,000  feet 
[1,  2,  3,  4,  5,  6].  Each  Center  is  sub-divided  into  about  20  sectors,  with  at  least  one  air  traffic 
controller  responsible  for  each  sector.  Each  sector  controller  may  talk  to  25-30  aircraft  at  a 
given  time  (the  maximum  allowed  number  of  aircraft  per  sector  depends  on  the  sector  itself). 
The  controller  is  in  charge  of  preventing  losses  of  separation  between  aircraft,  keeping  them 
separated  by  more  than  5  nautical  miles  horizontally,  and  2000  feet  vertically.  In  general, 
the  controller  has  access  to  the  aircraft’s  flight  plan  and  may  revise  the  altitude  and  provide 
temporary  heading  assignments,  amend  the  route,  speed,  or  profile,  in  order  to  attempt  to 
optimize  the  flow  and  to  keep  aircraft  separated ,  as  well  as  to  provide  weather  reports  and 
winds.  An  illustration  of  the  current  control  structure  is  presented  in  Figure  1. 

Existing  NAS  modeling  tools  have  functionality  which  spans  the  modeling  of  runway  and 
airport  capacity  and  operations,  through  airspace  operations  and  conflict  resolution  [7,  8],  to 
human  factors  and  man-machine  integration.  See  [9,  10]  for  detailed  surveys  of  NAS  modeling 
and  conflict  detection  and  resolution  methods.  A  recent  tool,  FACET  [11,  12],  represents  the 
first  accurate  NAS  simulation  tool,  with  the  additional  capability  of  a  “playback”  mode  using 
actual  traffic  flow  data.  Our  goal  is  to  develop  a  model  which  complements  existing  tools,  by 
providing  a  control  theoretic  component  which  models  the  influence  of  Air  Traffic  Control 
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(ATC).  While  the  additional  logic  required  to  model  the  actions  of  the  controller  does  not 
pose  a  significant  computational  problem  if  the  aircraft  density  in  the  airspace  is  low,  it 
becomes  an  issue  as  the  density  increases.  The  long  term  goal  of  increasing  capacity  as  well 
as  safety  in  the  NAS  cannot  be  achieved  without  an  in-depth  analysis  of  the  applied  control 
logic.  If  such  a  system  were  shown  to  model  the  current  airspace  with  sufficient  accuracy, 
then  a  wide  array  of  applications  would  become  feasible,  including  providing  additional 
support  for  ATC  in  predicting  delay. 

In  this  paper  we  present  a  model,  as  well  as  analytic  and  simulation  results,  of  the 
aircraft  and  controller  actions  within  a  sector  of  airspace.  Our  approach  is  Lagrangian , 
meaning  that  it  accounts  for  the  trajectories  of  the  aircraft  and  parameters  transported 
along  them,  for  example,  local  aircraft  density  (average  number  of  aircraft  per  surface  unit 
for  a  given  portion  of  airspace),  momentum  or  average  speed.  Most  of  the  existing  NAS 
models  use  this  framework,  see  for  example  [13,  14,  15,  11,  16,  17,  18].  However,  we  bridge 
our  results  with  Eulerian  approaches  such  as  [19].  Eulerian  approaches  are  control  volume 
based:  they  account  for  temporal  fluctuations  of  quantities  in  a  given  volume,  for  example, 
the  number  of  aircraft  present  in  a  portion  of  airspace.  The  connection  between  these  two 
approaches  is  made  through  a  new  concept  which  we  introduce:  sector  dynamic  capacity , 
also  linked  to  other  Lagrangian  studies  [14]  and  complexity  metrics  [20,  21],  We  use  our 
model  to  study  the  effect  of  aircraft  flow  density  requirements  at  sector  boundaries,  due,  for 
example,  to  miles- in-trail  requirements  at  airports  in  subsequent  sectors:  the  model  allows 
for  prediction  capability  as  to  how  the  current  system  might  react  to  imposed  flow  conditions. 
In  addition,  the  model  allows  for  the  testing  of  the  effectiveness  of  different  controller  policies 
in  minimizing  delays  in  posted  flight  plans. 

This  article  has  two  components:  airspace  modeling  and  analysis,  and  validation  and 
simulation.  In  the  first,  we  present  a  mathematical  model  for  a  controlled  sector,  based 
on  a  hybrid  system  model  for  each  aircraft  which  encodes  simple  aircraft  dynamics  under 
the  discrete  action  of  the  controller.  We  have  observed  that  the  set  of  commands  used  by 
controllers,  while  large,  is  finite,  and  consists  of  simple  actions  such  as:  “turn  to  heading  of 
x  degrees”,  “hold  current  heading”,  “fly  direct  to  jetway  y'\  “increase  speed  to  z  knots”.  We 
then  analyze  this  model  and  define  the  concept  of  sector  dynamic  capacity,  which  we  com¬ 
bine  with  analysis  to  predict  the  time  it  takes  to  overload  given  sectors  of  airspace,  and  thus 
predict  delays,  assuming  controllers  use  a  subset  of  their  available  control  actions.  We  finally 
link  these  results  with  Eulerian  approaches.  In  the  second  part  of  this  article,  we  validate 
the  previous  results  against  real  data.  Since  our  results  cannot  be  tested  on  the  real  ATC 
system  directly,  we  design  a  simulator  of  the  system,  which  we  implement  in  c-| — \-  interfaced 
with  MATLAB.  This  simulator  consists  of  the  mathematical  model  derived  in  the  first  part 
of  the  paper,  augmented  with  a  logic  for  switching  between  the  different  controller  states. 
This  logic  models  the  actions  of  a  human  air  traffic  controller,  and  results  from  minimizing 
a  particular  cost  function,  which  we  derive.  We  validate  this  simulator  against  Enhanced 
Traffic  Management  System  data  (ETMS)  and  then  show  that  our  analytical  predictions  for 
sector  capacity  are  effectively  observed  in  simulations.  Finally,  we  use  our  simulator  to  iden¬ 
tify  flow  conditions  under  which  conflict  resolution  decorrelates  from  metering  problems  (i.e. 
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scheduling  aircraft  according  prescribed  rules).  This  result  has  implications  for  numerous 
ATC  flow  management  studies  which  rely  implicitly  on  this  assumption  [22,  23,  19,  14]. 

The  data  that  we  present  in  this  paper  pertain  to  several  sectors  within  the  Oakland 
ARTCC,  located  in  Fremont,  CA.  Our  methodology,  however,  is  general  and  would  apply 
to  any  other  En- Route  portion  of  the  NAS.  The  geographical  data  presented  is  available  in 
the  public  domain  (we  used  JEPPESEN  [24]  high  altitude  en  route  charts).  Our  controller 
model  and  cost  function  have  been  designed  based  on  our  observations,  over  several  hours,  of 
sector  controllers  for  given  sectors;  we  have  made  justified  approximations  where  appropriate 
for  the  study  in  this  paper.  While  our  model  is  general,  most  of  the  scenarios  considered  in 
this  paper  do  not  represent  normal  traffic  flow:  because  we  are  interested  in  modeling  delay 
propagation  of  the  system  under  stress,  the  traffic  scenarios  modeled  represent  heavy  traffic 
flow  along  jetways. 

The  contributions  of  this  paper  are  a  new  mathematical  model  for  airspace  sectors,  which 
uses  hybrid  system  theory,  an  analytical  solution  to  the  Lagrangian  problem  of  delay  prop¬ 
agation  in  the  network,  and  its  link  with  Eulerian  approaches,  and  the  concept  and  use  of 
sector  dynamic  capacity.  On  the  application  side,  the  novelty  is  in  the  validation  of  a  control 
theoretic  model  of  the  human  air  traffic  controller,  as  well  as  in  the  validation  of  the  analyti¬ 
cal  predictions  against  real  data.  Finally,  the  decorrelation  results  shown  by  our  simulations 
are  new. 

This  paper  is  organized  as  follows:  Section  2  presents  the  model  used  for  aircraft  dynamics 
and  air  traffic  controller  actuation  (Section  2.1).  This  model  is  used  in  Section  2.2  to  predict 
the  propagation  of  airspace  congestion,  and  to  define  sector  capacity.  Section  3  presents  the 
design  of  the  simulator  (Section  3.1),  its  use  in  validating  our  analytical  predictions  (Section 
3.3),  and  its  use  in  demonstrating  the  decorrelation  between  conflict  resolution  and  flow 
metering  (Section  3.4). 

2  Air  Traffic  Flow  modeling  and  analysis 

The  structure  of  the  NAS  is  complex,  for  it  involves  a  multitude  of  interacting  agents 
and  technologies:  aircraft  monitoring,  flow  management,  communication,  and  human-in- 
the-loop.  For  the  present  work,  in  which  we  are  interested  in  predicting  delay,  we  extract 
and  model  only  the  features  which  are  important  for  this  purpose.  We  model  a  portion  of  the 
Oakland  ARTCC,  which  contains  five  sectors.  These  sectors  surround  the  Oakland  TRA- 
CON  (Terminal  Radar  Approach  CONtrol),  which  controls  the  aircraft  on  their  approaches 
into  San  Francisco,  San  Jose  and  Oakland  airports.  The  TRACON  is  the  final  destination 
of  the  traffic  that  we  consider  in  this  paper. 

We  model  a  sector  by  a  portion  of  airspace  containing  aircraft  under  the  local  control  of 
the  responsible  air  traffic  controller  (Figure  2).  The  interior  of  this  domain  is  the  controlled 
area  (in  which  the  local  controller  can  actuate  the  flow).  Within  each  sector,  navigation 
infrastructure,  including  jetways,  waypoints  and  navigation  aids ,  is  used  to  help  the  flow 
follow  desired  patterns;  we  therefore  include  it  in  the  model  and  use  it,  even  if  it  is  observed 
that  more  than  40%  of  the  aircraft  may  fly  off  the  jetway  at  any  given  time.  Our  model  allows 
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for  aircraft  to  fly  at  different  altitudes,  but  not  to  climb  or  descend.  Altitude  changes  are 
not  crucial  for  the  effects  we  want  to  identify:  the  type  of  sector  overloads  we  are  interested 
in  mainly  result  from  aircraft  acceptance  rates  at  destination  airports. 


2.1  Aircraft  behavior 


We  use  a  hybrid  model  for  each  aircraft.  A  hybrid  model  describes  the  evolution  of  a  system 
by  a  set  of  discrete  modes ,  each  associated  with  a  continuous  dynamical  system ,  and  discrete 
switches  which  enable  the  system  to  jump  from  one  mode  to  another  instantaneously.  In 
mathematical  terms,  we  will  describe  the  motion  of  aircraft  i  by: 


x 


i 


dxi 

dt 


^current  heading 


(1) 


where  ucurrent  heading  G  R2  is  a  constant  velocity  vector  held  by  the  aircraft  until  the  next 
discrete  switch:  a  heading  or  speed  change  which  changes  ^current  heading-  G  R2  is  the 
planar  position  of  aircraft  i.  Integration  of  equation  (1)  over  time  produces  a  continuous 
piecewise  affine  trajectory.  We  prefer  such  a  model  over  a  continuous  dynamical  model  for 
two  reasons.  First,  the  time  scale  of  a  “change  in  aircraft  behavior”  (for  example  a  turn  or 
slow  down)  is  on  the  order  of  30  seconds,  whereas  the  time  scale  of  a  straight  line  portion 
of  the  flight  is  usually  much  longer,  sometimes  half  an  hour  or  more,  thus  we  ignore  the 
dynamics  of  such  maneuvers  and  focus  on  their  effects  only  (the  set  of  resulting  straight 
lines).  Second,  the  update  rate  of  ATC  monitoring  is  in  general  not  more  than  30  seconds, 
which  makes  the  details  of  these  maneuvers  inaccessible  to  the  ATC.  This  approximation  is 
widely  accepted  in  literature  [25,  26,  27,  17,  28,  22], 

Monitoring  ATC  shows  that  a  finite  set  of  maneuvers,  which  depends  on  local  parameters, 
is  used.  Combinations  of  these  maneuvers  result  in  a  conflict-free  flight  environment  in  which 
the  constraints  of  the  air  traffic  flow  are  met.  The  maneuvers  shown  in  Figure  4  are  modeled 
(and  consist  of  changing  the  right  hand  side  of  (1)  according  to  certain  rules  which  we  now 
make  explicit).  The  validity  of  models  similar  to  this  has  been  confirmed  by  statistical 
studies  realized  at  MIT  ICAT  (Histon  and  Hansman  [20]). 


1.  Speed  change:  ATC  may  decelerate  or  accelerate  the  aircraft  along  its  flight  plan: 

^modified  speed  •  A  ^current  heading  (2) 

where  A  G  R+  defines  the  magnitude  of  the  velocity  change.  Our  model  will  allow  a 
finite  set  of  speeds  (which  means  A  has  a  finite  number  of  acceptable  values).  This 
encodes  the  fact  that  the  ATC  generally  has  a  finite  set  of  possibilities  in  the  choice  of 
speeds,  because  the  aircraft  flies  at  its  optimal  speed  per  altitude  and  ATC  will  speed 
up  or  slow  down  the  aircraft  by  not  more  than  10%  of  the  current  value. 

2.  Vector-for-spacing  (VFS):  This  maneuver  consists  of  a  deviation  of  the  aircraft  from  its 
original  flight  plan  for  a  short  time  (part  1  of  the  maneuver),  and  a  second  deviation, 
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bringing  it  back  to  its  original  flight  plan  (part  2  of  the  maneuver).  This  stretches  the 
path  that  the  aircraft  must  follow,  and  therefore  generates  a  delay.  The  length  of  this 
maneuver  depends  on  the  geometry  of  the  sector.  Calling  the  rotation  matrix  by 
angle  ip,  we  have: 

Apart  1  :=  Rp  ■  ^current  heading  (First  half  of  the  maneuver)  , 

Apart  2  :=  R- 24>  •  Fpart  i  (Second  half  of  the  maneuver)  ' 

3.  Shortcut  /  Detour:  In  certain  situations,  the  ATC  will  have  the  aircraft  “cut”  between 
two  jetways,  a  maneuver  which  could  either  shorten  or  lengthen  the  flight  plan.  The 
decision  to  command  such  a  maneuver  is  often  dictated  by  conflict  resolution,  but  could 
also  be  used  to  shorten  the  overall  flight  time  if  sector  occupancy  allows  it  (sometimes 
called  “direct-to”  by  pilots): 

^shortcut  :=  R-0  ■  ^current  heading  for  the  duration  of  the  maneuver  (4) 

until  the  next  ATC  action  is  taken.  Here  again  ip  is  the  angle  by  which  ATC  turns  the 
aircraft  to  achieve  the  shortcut. 

4.  Holding  pattern:  In  some  extreme  conditions,  holding  patterns  are  used  to  maintain  an 
aircraft  in  a  given  region  of  space  before  eventually  letting  it  follow  its  original  flight 
plan.  This  is  modeled  by  assigning  the  aircraft  to  a  predefined  zone  and  keeping  it 
there  while  preventing  other  aircraft  from  entering  that  zone. 

2.2  Lagrangian  analysis  of  delay  propagation  in  the  NAS 

A  large  proportion  of  air  traffic  jams,  i.e.  portions  of  airspace  saturated  by  aircraft,  is  gen¬ 
erated  by  restrictions  imposed  at  destination  airports,  usually  themselves  driven  by  weather 
or  airport  congestion.  These  restrictions  are  imposed  as  either  miles-in-trail  or  minutes-in- 
trail,  representing  the  distance  (or  time)  required  between  aircraft  in  a  flow  incoming  to  the 
TRACON,  and  are  referred  to  as  metering  constraints.  Figure  5  illustrates  the  topology  of 
the  flows  incoming  to  San  Francisco  Airport  (SFO)  which  are  often  subject  to  this  type  of 
constraint.  These  constraints  tend  to  propagate  backwards  from  the  airport  into  the  net¬ 
work,  and  result  in  miles-in-trail  constraints  imposed  at  the  entry  points  of  each  sector.  For 
example,  in  the  case  of  Figure  5,  typically,  the  backpropagation  of  these  metering  conditions 
is  as  follows:  TRACON  — >  sector  34  — >  sector  33  — >  Salt  Lake  Center  •  •  •  and  similarly  for 
the  two  other  flows. 

In  the  current  system,  these  restrictions  are  imposed  empirically.  We  would  like  to 
understand  (%)  how  the  traffic  jams  propagate;  and  (ii)  what  the  optimal  control  policy  should 
be  under  these  restrictions,  in  order  to  ensure  maximal  throughput  into  the  TRACON. 

2.2.1  Shock  wave  propagation 

We  present  a  simple  Lagrangian  model  of  merging  flows,  introduced  earlier  in  [29,  30].  This 
model  of  merging  flows  predicts  the  backpropagation  of  a  traffic  jam  from  a  destination 
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airport  into  the  network.  It  is  Lagrangian,  because  it  models  the  trajectories  of  all  agents 
in  the  considered  space,  rather  than  averaged  quantities  in  a  control  volume,  such  as  the 
number  of  aircraft  in  a  given  sector.  This  model  can  be  applied  to  the  merging  flows  of 
the  type  shown  in  Figure  5.  We  will  use  it  here  to  derive  the  dynamic  capacity  of  a  sector. 
Given  prescribed  inflow  and  outflow  conditions,  we  define  the  dynamic  capacity  of  a  sector 
to  be  the  number  of  aircraft  which  can  be  actuated1  in  this  sector,  until  ATC  action  has  to 
occur  upstream  from  this  sector  in  order  to  not  violate  the  outflow  conditions.  We  assume 
in  this  definition  that  the  sector  is  initially  empty,  but  we  will  show  how  the  definition  may 
be  applied  to  cases  in  which  the  sector  is  initially  non-empty.  The  dynamic  capacity  is  a 
concept  which  appears  naturally  in  the  following  problem: 

Given  a  required  spacing  between  the  aircraft  ( metering  constraint)  of  ATout,  compute  a 
controller  policy  which  will  force  groups  of  aircraft  to  exactly  satisfy  the  metering  constraint 
at  the  sector  exit  point  (each  aircraft  separated  by  exactly  ATout)  while  maintaining  separation 
at  all  times. 

We  first  consider  a  very  simple  version  of  the  problem,  in  which  the  controller  uses 
only  two  modes  (fast  and  slow).  This  model  is  not  overly  restrictive:  we  can  use  several 
methodologies  to  map  the  full  automaton  of  Figure  4  to  this  model:  see  for  example  [23,  14]. 

We  introduce  the  initial  arc  length  distance  a)  e  R  of  aircraft  i  along  its  arrival  route  to 
the  airport.  For  example,  a)  =  —200  means  that  aircraft  i  has  to  fly  200  nm  before  landing 
in  the  airport.  We  call  xex  G  R  the  location  at  which  the  metering  condition  is  imposed. 
For  example,  xex  =  —50  means  that  the  metering  is  applied  50  nm  from  the  airport.  It 
is  possible  to  assume  without  loss  of  generality  that  the  aircraft  are  numbered  in  order  of 
arrival  (the  af)  are  indexed  in  increasing  order  of  magnitude). 

We  consider  the  following  situation:  all  aircraft  are  initially  at  maximum  speed  vmax , 
and  in  order  to  enforce  metering,  ATC  slows  down  aircraft  i  to  its  minimum  speed  nmin  (see 
Figure  6),  at  a  location  a;fvltch  and  time  f®wltch  which  we  will  determine  (see  Figure  7).  This 
scenario  is  represented  as  a  dash-dot  line  in  Figure  4.  We  impose  the  condition  that  each 
aircraft  cross  the  metering  point  xex  at  exactly  fbiock  +  (i  —  l)ATout  where  fbiock  is  the  time 
at  which  the  metering  condition  is  initiated.  For  simplicity,  we  take  the  origin  of  time  at 
t0  =  0.  This  leads  to  the  following  kinematic  equations  of  the  aircraft  under  this  actuation: 

Xi(t)  =  a°  +  vmaxt  t  e  [0,  f -witch] 

Xi(t)  =  hi  +  vminf  t  E  [f™tch,  fbiock  +  (*  -  1)  ATout] 

The  assumption  of  continuity  of  ay(f)  enables  us  to  solve  for  bt  J  from  which  we  deduce  the 
switching  time  as:  f®wltch  =  (bj  —  a ))/ (vmxr  —  iu).  Under  the  following  feasibility  conditions: 

^  [*^ex  ^max  (tbloc.k  (f  1)  ATout) ,  Xex  Umin(tblock  if  l)ATout;)]  (5) 

the  propagation  speed  of  the  traffic  jam  may  be  computed  analytically,  by  solving  for  the 

1ATC  actuation  of  an  aircraft  means  any  alteration  of  its  original  flight  plan,  in  order  to  meet  desired 
conditions  in  the  region  of  interest. 


location  of  the  edge  of  the  traffic  jam  in  space  and  time: 


^.switch  _ 

Li 

^  switch  _ 


%ex  'Rmin^block  1)AZ/  CL 9 

h'max  h?min 

0  |  %ax[^ex_ ^min^block- (^—  1)AL  — 0^  ] 
Cl  „■  — | 

1  h’rnax- h?min 


(6) 


where  AT  :=  uminATout  is  the  metered  spacing  at  the  outflow  of  the  sector.  At  a  given  time 
t,  we  call  traffic  jam  or  metered  platoon  the  set  of  aircraft  such  that  t®wltch  <  t.  These  aircraft 
have  already  been  actuated  (see  Figure  6).  It  follows  directly  from  (6)  that  the  traffic  jam 
will  not  grow  if  the  two  following  conditions  are  met: 


/switch  ^  -/switch 
Li  "  Li+ 1 

...switch  ^  switch 
xi+ 1  ^  xi 


AL  <  a°  -  a°+1 


(7) 


Condition  (6)  above  is  a  sufficient  condition  for  the  length  of  the  traffic  jam  to  decay  (which 
can  be  observed  by  inspection  of  the  slope  of  the  switching  curve  or  shock  wave  of  points 
(xfwltch ,f®wltch )  displayed  in  Figures  7  and  8).  It  is  a  local  property  of  the  problem:  the 
conditions  depend  only  on  af  —  a°+1,  not  on  all  aircraft  considered  here.  The  second  equation 
in  (6)  is  in  fact  a  one-dimensional  discretized  steady  Lighthill-Whitham-Richard  (LWR) 
equation,  which  appears  naturally  in  highway  congestion  problems  [31].  This  fact  links  this 
Lagrangian  approach,  based  on  aircraft  trajectory  analysis,  with  Eulerian  approaches  such 
as  [19],  which  are  based  on  conservation  equations;  it  relates  local  properties  of  the  flow 
(local  monotonicity  of  a  variable,  here  the  ^-location  of  the  wavefront)  to  global  quantities 
(here  the  trajectories  of  the  aircraft).  This  resiilt  is  illustrated  in  Figure  7. 

The  condition  that  each  aircraft  reaches  xex  exactly  at  the  time  prescribed  is  restrictive, 
what  is  important  is  the  flow  rate  but  not  the  actual  crossing  times.  Therefore  it  would  be 
of  greater  use  to  pose  the  problem  as  follows:  Given  {a/)}ie[i.jvp  compute  the  switching  policy 
which  delivers  at  most  one  aircraft  every  A Tout  sec.  at  the  location  xex  while  maintaining 
separation,  and  minimizes  the  arrival  time  of  aircraft  N .  This  problem  may  be  posed  as 
a  linear  program:  (a.)  minimize  the  arrival  time  of  aircraft  N,  (b.)  while  separating  the 
aircraft  by  more  than  ATollt  at  xex,  (c.)  with  at  most  one  switch  between  the  initial  position 
a°  <  xex  and  the  exit  xex  of  the  considered  airspace: 


a.  Minimize  [0,  •  •  ■  0,  — 1]6 

b.  Subject  to 


■  -1 

1 

0 

...  o' 

A  Tout 

0  - 

-1 

1 

-1 

1  0 

b  CZ  ^min 

0  • 

0 

-1  1 

A  Tout 

akH 

^min 

if 

+  b 

„  .  > 
^min 

,1]T 

Umax 

■Rmax  ) 

where  a  =  [af  •  •  •  ,  a°N]T  and  b  =  [6i,  •  *  -  ,  &at]t.  Note  that  the  right  hand  side  of  (b.)  can  be 
changed  to  [AXj,  •  •  •  ,  A7/v]T  in  order  to  account  for  time- varying  metering  conditions.  The 
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advantages  of  this  formulation  are  that  one  includes  the  possibility  to  optimize  an  objective 
function,  which  is  in  this  case  the  arrival  time  of  the  last  aircraft  in  the  platoon,  and  that 
one  is  able  to  deal  with  time  varying  acceptance  rates  at  the  airport.  Conditions  (7)  for 
shock  monotonicity  derived  previously  are  still  valid  locally  for  any  solution  derived  with  the 
linear  program  above. 


2.2.2  Sector  overload  predictions 

Using  the  analysis  of  the  previous  section,  it  is  fairly  easy  to  predict  the  dynamic  capacity 
of  a  sector.  Consider  the  worst  case  scenario:  an  incoming  flow  of  aircraft,  each  at  umax , 
separated  in  time  by  AT;n  chosen  to  violate  the  second  condition  in  (7).  This  will  create  a 
traffic  jam  originating  at  SFO  (equated  with  entrance  to  SFO  TRACON),  which  “piles  up” 
and  progressively  fills  sector  34.  We  call  l  the  arc  length  distance  of  the  portion  of  the  arrival 
jetway  contained  in  sector  34  (see  Figure  6).  We  assume  that  the  sector  is  initially  empty. 
Using  equations  (6),  we  can  compute  the  maximum  number  iViimit  of  aircraft  which  can  be 
stacked  along  the  length  l  of  the  jetway  in  the  sector,  before  no  space  is  available  anymore 
on  this  jetway.  These  aircraft  are  labeled  as  metered  platoon  in  Figure  6;  for  example,  for  the 
situation  shown  in  this  Figure,  approximately  half  of  l  is  occupied  by  the  metered  platoon 
at  the  time  considered,  so  the  number  of  metered  aircraft  is  approximately  Wimit/2.  When 
the  number  of  metered  aircraft  reaches  Wimit,  after  a  time  called  Tiimit. ,  the  actuation  shown 
in  Figure  6  (ATC  slows  aircraft  down)  has  to  occur  upstream  (in  sector  33).  Wimit  and  Tiimit 
are  respectively  given  by: 


I  (Unax  Unin) 


Unax  Unin  (  ^7>ut  ATjn ) 


(8) 


y  _ 

limit 


^max^-^in 


out 


7max  ^min 


AT 


out 


AT,, 
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If  aircraft  are  initially  present  in  the  sector,  these  two  quantities  can  be  modified  by  replacing 
l  by  the  distance  of  the  last  aircraft  arrived  in  the  sector  to  the  entrance  of  the  sector.  We 
call  IViimit  dynamic  capacity  (and  not  static)  because  it  depends  on  the  inflow  and  outflow 
conditions  (and  not  only  on  geometric  parameters). 

Several  comments  can  be  made  regarding  the  two  previous  results:  (%)  As  Unin  —  Unax  — >  0, 
IViimit.  — >  0:  no  aircraft  can  be  handled  in  the  sector  because  no  actuation  is  possible  (it  is  not 
possible  to  “make  the  aircraft  lose  time”  in  this  sector);  (ii)  As  ATout  — ATin  — >  0,  IViimit.  — >  oo 
and  Tiimit  — >  oo:  if  the  incoming  flow  is  such  that  it  is  almost  metered  as  imposed  at  the  exit 
of  the  sector,  the  number  of  aircraft  required  to  saturate  this  airspace  becomes  large  and  the 
time  it  takes  to  saturate  this  sector  grows  accordingly. 

The  construction  of  the  switching  curve  (x?wltch;  ^?wltch)  described  previously,  can  be  used 
to  compute  the  maximal  extent  of  a  traffic  jam,  or  the  portion  of  jetway  affected  by  a  traffic 
jam.  Using  (6),  one  can  trace  the  shock  location  in  the  (x,t)  plane  (Figure  8).  The  edge  of 
the  traffic  jam,  called  xm,  obtained  at  tm  gives  the  worst  situation  obtained  from  the  initial 


10 


configuration  a®  of  the  aircraft:  at  tm,  the  extent  of  the  traffic  jam  is  at  its  maximum.  In  the 
case  of  Figure  8,  we  see  that  the  traffic  jam  does  not  propagate  more  than  300  nm  upstream 
from  the  destination  of  the  aircraft,  called  xex.  Therefore  no  metering  conditions  need  to  be 
applied  upstream  from  that  point.  In  the  current  system,  such  information  is  not  available 
to  the  ATC,  thus  leading  to  extra  buffers  taken  by  the  controllers,  which  in  turn  leads  to 
non  optimal  operating  conditions  as  well  as  backpropagation  of  “virtual  overload” ,  a  set  of 
conservative  precautions. 

3  Validation  against  ETMS  data 

The  models  of  the  previous  section  rely  on  a  mathematical  analysis  of  metering,  based  on 
geometry.  In  order  to  understand  how  realistic  these  models  are,  we  need  to  validate  them 
against  the  real  system.  Since  the  real  system  is  not  available  as  a  testbed,  we  have  designed 
an  abstraction  of  it  which  reproduces  its  behavior  adequately.  In  this  section,  we  present  the 
design  of  a  simulator  which  mimics  true  ATC  behavior,  and  the  validation  of  this  simulator 
against  ETMS  data.  The  graphhical  interface  of  the  simulator  is  shown  in  Figure  3.  In 
Section  3.3,  we  validate  the  analytical  predictions  of  Section  2,  using  the  simulator,  and  in 
Section  3.4,  we  show  how  the  simulator  may  be  used  to  derive  conditions  on  the  decorrelation 
of  flow  metering  and  conflict  resolution. 

3.1  Simulator  design 

Our  simulator  is  based  on  empirical  studies  that  we  realized  at  the  Oakland  ARTCC,  its  core 
is  based  on  observed  behavior.  Figure  4  (all  transitions  enabled)  summarizes  the  behavior 
model  observed  at  the  ARTCC.  The  switching  logic  behind  the  transitions  is  the  object  of 
this  section,  and  is  implemented  in  the  form  of  a  cost  function,  which  we  describe  in  this 
section. 

3.1.1  Overall  Program  Flow 

The  overall  program  flow  of  the  simulator  is  shown  in  Figure  9.  The  input  to  the  code  is  a 
set  of  aircraft  filed  flight  plans  (Figure  9,  middle  column),  that  can  either  be  user  generated 
or  taken  from  ETMS  data  (as  in  FACET).  As  in  the  true  system,  these  flight  plans  are  not 
conflict-free  and  usually  do  not  satisfy  metering  conditions  imposed  on  the  network.  Once 
the  program  is  initialized,  aircraft  motion  simulation  follows  these  flight  plans  (Figure  9,  left 
column).  As  time  is  advanced,  conflict  as  well  as  metering  constraints  are  dealt  with  on  a 
sector  by  sector  basis  (with  sector-wide  look  ahead,  Figure  9,  right  column),  according  to 
the  full  automaton  shown  in  Figure  4.  The  flight  plans  are  updated  accordingly. 
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3.1.2  Key  Data  Structures 


Aircraft  dynamic  equations  (1)  produce  a  set  of  segments;  the  knowledge  of  the  points 
connecting  the  segments  and  of  the  aircraft  velocity  is  thus  enough  to  define  an  aircraft 
trajectory.  This  trajectory  is  thus  implemented  as  a  linked  list  of  points  [x,y,z],  with  a 
prescribed  velocity  between  the  points.  The  linked  list  is  modified  by  the  simulated  controller 
in  the  program.  The  output  for  each  aircraft  is  the  updated  linked  list.  The  sectors  are 
implemented  as  sets  of  linked  lists  accessible  to  a  controller.  They  also  contain  additional 
data  such  as  metering  conditions  (number  of  aircraft  through  a  given  boundary  per  time 
unit). 


3.1.3  Controller  Emulation 

ATC  behavior  is  modeled  by  three  levels  of  priority: 

•  Priority  1:  No  loss  of  separation.  The  prevalent  requirement  for  ATC  is  to  ensure  that 
any  aircraft  pair  is  always  separated  by  more  than  5  nautical  miles. 

•  Priority  2:  Metering  conditions.  The  controller  needs  to  ensure  that  the  outflow  from 
his  sector  is  an  acceptable  inflow  for  the  next  sector  (or  TRACON).  Metering  conditions 
can  be  of  various  nature:  admittance  rate  or  separation  at  downstream  junctions. 

•  Priority  3:  Best  possible  throughput.  ATC  will  try  if  possible  to  give  direct  routes  to 
aircraft  in  order  to  minimize  their  flight  times. 

These  priorities  may  be  modeled  using  the  cost  function  J\ 


J 


costLOs  +  costBc  breach  T  COStjgiay  T  COSta;rcraft  actuation 


+  cost 


maneuver 


+  cost 
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(10) 


Each  term  of  the  cost  is  a  weighted  function: 


^  =  E 
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1.  Loss  of  separation  (LOS)  cost:  nlLOS  is  the  number  of  losses  of  separation  involving 
aircraft  i  in  the  current  sector  with  its  current  flight  plan.  A Tfos  is  the  time  until  the 
first  loss  of  separation  for  aircraft  i. 

2.  Boundary  condition  (BC)  breach  cost:  T1(reach  is  the  time  by  which  an  aircraft  violates 
the  AT  time  separation  constraint  from  its  predecessor  (set  to  zero  if  the  two  aircraft 
are  separated  by  more  than  AT). 
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3.  Delay  cost:  TOA^-TOA^  accounts  for  the  difference  between  predicted  and  actual 
time  of  arrival  (TO A)  at  the  last  waypoint  of  the  flight.  Positive  delays  are  penalized; 
earlier  arrivals  are  favored.  TOAlvred  and  TOA\eai  are  computed  by  integration  of  the 
flight  plans  for  each  aircraft. 

4.  Aircraft  actuation  cost:  Nmoved  accounts  for  the  number  of  flight  plan  modifications 
chosen  in  the  current  solution.  Large  Nmoved  are  penalized  (the  solution  chosen  by  the 
ATC  is  often  the  simplest). 

5.  Maneuver  cost:  Maneuver  accounts  for  the  cost  of  the  maneuver  selected  for  aircraft  i. 

Not  all  maneuvers  are  of  equal  preference  and  therefore  have  different  costs.  It  is  as 
easy  for  a  controller  to  prescribe  a  10%  speed  change,  a  VFS  or  a  shortcut.  A  holding 
pattern  is  the  least  preferred  option,  for  it  requires  constant  monitoring  of  the  aircraft. 
This  reflects  in  the  weight  choice:  Jlspeed  change  ~  Shortcut  ~  4fS  <  folding  pattern' 

The  ratio  folding  pattern/' 4>eed  change  is  °n  the  °rder  °f  10' 

6.  Minimal  distance  cost:  f(dlmin )  penalizes  aircraft  distributions  in  which  aircraft  are 
closely  spaced  (but  do  not  lose  separation)  against  more  sparse  distributions.  Here, 
distmax  =  7nm. 

/(4nn)  =  dh  ■  Wdist  if  < in  <  clist  max 

inin 

f(d  ^n)  =  0  otherwise. 

In  order  to  reflect  the  three  levels  of  priority  of  the  human  ATC  stated  earlier,  the  weights 
shown  in  the  cost  function  J  are:  Wlos  rs_/  lO300  »  ^breach  104  S>  other  weights  ~  10. 
Thus  a  computation  which  tries  to  minimize  J  will  hrst  deal  with  losses  of  separation,  then 
metering  conditions,  and  finally  optimization  of  the  flow.  An  example  of  a  cost  landscape 
for  a  given  topology  of  two  aircraft  is  illustrated  in  Figure  10. 

In  order  to  reduce  the  computational  time,  we  define  as  AAoice  the  maximum  number 
of  aircraft  considered  for  maneuvering  by  the  simulated  controller  in  each  time  iteration. 
We  restrict  A^hoice  to  be  generally  less  than  the  actual  number  of  aircraft  per  sector.  This 
term  is  a  trade-off  between  run-time  and  control-quality  and  in  our  simulations  it  was  set  in 
the  range  of  4  to  8,  depending  on  the  targeted  goals.  Aircraft  are  selected  according  to  the 
following  rule:  aircraft  involved  in  LOS  are  selected  hrst,  then  aircraft  breaching  boundary 
conditions,  and  finally  remaining  aircraft  until  the  selection  list  has  reached  Nc h0ice  aircraft, 
or  until  there  are  no  more  aircraft  to  select.  In  practice,  4  <  AAoice  <  8,  where  A^hoice  =  8 
enables  more  complicated  situations  but  makes  the  code  run  more  slowly.  The  set  of  all 
maneuver  combinations  for  the  AT^oice  aircraft  is  called  the  maneuver  set. 

At  each  iteration  of  the  controller  emulation  loop,  an  exhaustive  recursive  search  on  the 
maneuver  set  is  run  for  the  chosen  aircraft  in  order  to  find  a  set  of  AAioice  maneuvers  which 
minimizes  J.  The  computational  complexity  of  finding  the  optimal  J  for  A^hoice  aircraft 
subject  to  nmaneuver  possible  discrete  maneuvers  is  0((nmaneUver)Archoice)-  We  can  reduce  this 
cost  to  O ( (nmaneUver  —  2)iVchoice):  (i)  the  cost  of  the  current  maneuver  has  already  been 
computed  at  the  previous  step  and  thus  does  not  need  to  be  recomputed;  (ii)  two  maneuvers 
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are  mutually  exclusive  (shortcuts),  therefore  only  one  needs  to  be  called.  Including  the  cost 
of  checking  for  conflicts,  the  total  cost  of  a  time  iteration  becomes: 

O  { N2  ■  (n  —  o)Nchoice) 

'S  ^iVmax  V  ‘'maneuver  ^  )  J 

where  7Vmax  represents  the  total  number  of  aircraft  in  the  sector.  Due  to  both  the  discretiza¬ 
tion  of  time,  and  the  restriction  of  the  search  space  to  a  manageable  number  of  aircraft, 
our  search  is  not  guaranteed  to  find  the  global  optimum.  However,  as  we  show  in  the  next 
section,  the  search  does  provide  a  reasonable  approximation  of  the  controller’s  behavior.  By 
adjusting  the  two  key  control  parameters,  the  number  of  selected  aircraft  Nc hoice  and  time 
between  controller  activation  A Tact,  a  transparent  trade-off  between  run-time  and  control 
quality  can  be  made. 

3.2  Code  validation  against  ETMS  data 

The  controller  logic  presented  in  the  previous  section  is  the  result  of  numerous  observations 
we  made  in  the  Oakland  ARTCC,  monitoring  the  work  of  air  traffic  controllers.  The  fact  that 
one  can  classify  ATC  action  into  a  set  of  preferred  directives  was  experimentally  validated 
for  a  different  airspace  (see  Hansman  and  Histon[20]).  However,  even  if  the  automaton  of 
Figure  4  and  the  cost  function  of  the  previous  section  implemented  in  the  simulator  are 
consistent  with  our  observations,  there  is  no  a  priori  guarantee  that  these  would  produce 
the  same  effects  on  the  system  as  a  human  controller.  For  this  reason,  we  assess  how  well 
our  cost  function  describes  the  decision  making  of  a  human  controller  by  validating  our  code 
against  recorded  aircraft  trajectories.  We  use  ETMS  data  provided  by  NASA  Ames.2 

We  first  explain  the  data  extraction  process  which  enables  us  to  convert  ETMS  data  to 
a  readable  format  for  our  simulator.  We  then  present  the  three  types  of  validations  realized. 

3.2.1  Data  extraction 

We  can  extract  two  types  of  information  from  the  ETMS  data:  the  actual  flown  aircraft 
trajectories,  and  the  hied  flight  plans  for  each  aircraft  (eventually  updated  if  modifications 
are  made  during  the  flight).  The  position  is  given  in  latitude  /  longitude  and  in  terms  of 
navaids,  fixes,  and  jetways  which  we  represent  in  cartesian  coordinates,  an  approximation 
valid  for  the  portion  of  airspace  of  interest  to  us.  The  hied  flight  plan  is  given  in  terms 
of  navaids,  hxes,  jetways,  which  we  also  convert  into  cartesian  coordinates  using  a  public 
database  (http://www.airnav.com).  Future  versions  of  our  simulator  might  use  recently 
developed  ETMS  analyses  tools  such  as  [13]  to  perform  these  tasks  automatically.  This  study 

2Data  is  collected  from  the  entire  population  of  flights  with  filed  flight  plans  in  the  NAS.  ETMS  data 
is  sent  from  the  Volpe  National  Transportation  System  Center  (VNTSC)  to  registered  participants  via  the 
Aircraft  Situation  Display  to  Industry  (ASDI)  electronic  file  server.  A  file  containing  all  recorded  data  is 
generated.  It  displays  for  each  aircraft  the  current  flight  data  (time,  position,  speed,  heading),  as  well  as 
the  filed  flight  plan  (in  terms  of  navaids,  jetways,  fixes,  etc.).  The  update  rate  of  the  measurements  is  of  the 
order  of  one  minute. 
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is  limited  to  sector  control,  and  we  did  not  implement  the  Traffic  Management  Unit  (TMU) 
action.  TMU  operates  at  the  Center  level  and  makes  strategic  flow  scheduling  decisions, 
which  go  beyond  the  range  of  a  single  sector  controller.  We  therefore  need  to  validate  our 
simulator  at  a  scale  at  which  TMU  actuation  is  already  incorporated  in  the  flight  plans 
(typically  one  or  two  sectors).  Since  our  interest  focuses  on  sectors  32,  33,  34,  15,  and  13, 
the  actual  flight  plans  are  truncated,  and  we  keep  only  points  in  that  sector.  The  filed 
flight  plans  are  truncated  similarly.  Their  estimated  time  of  entrance  in  the  sector  is  set  to 
the  actual  time  of  entrance  of  the  closest  actual  recorded  point.  The  entrance  point  in  the 
computational  domain  is  taken  to  be  the  closest  conflict-free  point  to  the  intersection  of  the 
flight  plan  and  the  boundary  of  the  corresponding  sector.  The  altitude  assigned  to  the  filed 
flight  plans  is  set  to  the  average  altitude  of  the  actual  trajectory  in  that  sector. 

3.2.2  Validation 

Comparison  of  flight  times.  We  select  the  first  100  aircraft  of  the  ETMS  data  set  which 
are  flying  above  33000ft  for  more  than  6  minutes.  Their  recorded  trajectories  are  extracted 
as  sequences  of  waypoints  which  are  used  as  initialized  flight  plans  for  our  simulations.  We 
then  compare  the  flight  times  in  the  simulation  to  the  actual  ones.  The  experiment  is  run 
for  the  following  set  of  speeds:  M  e  {0.6,  0.7,  0.8}  (M  is  the  Mach  Number),  which  matches 
our  observations  in  the  data  for  this  altitude.  In  the  run,  the  controller  is  activated  every 
ATact  =  lOsec.  The  results  are  shown  in  Figure  11.  Two  main  conclusions  can  be  made:  (i) 
the  simulator  is  able  to  recreate  the  flow  without  major  modification,  and  eventually  resolves 
apparent  conflicts  in  the  data  -  these  conflicts  can  be  due  to  inaccuracy  of  the  measurements 
or  transmission  (two  aircraft  separated  by  less  than  5  nm  at  the  same  altitude),  or  due  to 
problems  of  interpolation  when  speed  changes  in  time;  (ii)  the  time  comparison  (Figure  11) 
shows  relatively  good  matching.  The  flight  times  provided  by  the  simulator  are  usually 
shorter  than  in  reality  because  by  default  the  simulator  will  always  try  to  maximize  the 
throughput  in  the  sector.  The  mean  deviation  is  120  sec.  for  flight  plans  with  an  average 
duration  of  1300  sec.  (9.2%  error). 

Validation  of  conflict  resolution.  We  select  aircraft  flying  through  sector  33  in  a  time 
frame  of  10  hours  (a  total  set  of  314  aircraft)  and  simulate  these  flights.  The  filed  flight 
plans  are  not  conflict-free.  We  would  like  to  show  that  the  simulator  is  able  to  provide  a 
conflict-free  environment.  For  this  run,  the  activation  time  of  the  controller  is  A Tact  =  20sec.3 
The  set  of  speeds  allowed  is  M  e  {0.55,  0.75,  0.89}  (since  we  are  considering  the  full  range 
of  altitude,  we  need  to  consider  the  full  range  of  speed).  The  simulator  is  able  to  provide 
a  conflict  free  environment.  During  the  simulation,  it  has  to  actuate  50  different  aircraft. 
The  number  of  resolved  conflicts  can  be  assumed  to  be  on  the  same  order,  since  a  single 
intervention  will  usually  resolve  not  more  than  one  conflict. 

Validation  of  maneuver  assignments.  The  core  of  the  simulator  is  the  model  of  the  human 

3 A Tact  =  lOsec.  or  A Tact  =  20sec.  is  on  the  order  of  the  maximal  actuation  rate  of  a  controller.  We 
choose  ATact  =  20sec.  in  this  particular  case  because  of  the  duration  of  the  computation  (10  hours  of  real 
time  are  simulated). 
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controller  by  a  decision  procedure  based  on  the  cost  function  described  previously.  The 
validation  so  far  has  shown  the  correlation  of  flow  patterns  generated  by  our  code  and  these 
observed  in  reality.  We  would  also  like  to  assess  the  validity  of  our  decision  procedure.  For 
this,  we  identify  conflict  resolution  maneuvers,  which  are  typically  obtained  by  identifying 
deviations  from  the  hied  flight  plan.  For  these  maneuvers,  we  generate  the  following  input 
data:  all  aircraft  are  assigned  their  actual  flight  plan  (recorded  trajectories),  and  the  aircraft 
for  which  the  maneuver  is  identified  is  assigned  its  hied  flight  plan  (a  set  of  way  points). 
We  thus  put  the  simulator  in  the  same  situation  as  the  human  controller,  in  which  it  has 
to  make  the  decision  that  was  actually  taken.  For  sector  33,  we  were  able  to  identify  20 
distinct  maneuvers  out  of  300  examined  flight  plans.  The  simulator  reproduced  16  of  them. 
An  example  is  shown  in  Figure  12. 

3.3  Validation  of  the  analytical  predictions 

In  this  section,  we  assess  how  accurate  our  analytical  predictions  are  with  respect  to  the  real 
system,  by  comparing  them  with  simulations.  We  present  an  example  of  two  backpropagating 
shocks,  solved  with  an  extension  of  the  method  explained  in  Section  2.  Two  platoons,  each  at 
10  miles- in-trail5,  are  respectively  subjected  to  15  miles-in-trail  and  20  miles- in-trail  outflow 
boundary  conditions  (the  boundary  conditions  for  platoon  2  start  after  all  aircraft  of  platoon 
1  have  reached  the  TRACON  at  time  t  =  4300).  The  speeds  are  M  e  {0.59,  0.75,  0.89}. 
The  resulting  aircraft  flows  for  this  analytic  solution  are  shown  in  Figure  13.  The  results 
and  interpretations  are  shown  in  Figures  15  and  14.  Two  shocks  appear  successively.  From 
Figure  15,  we  can  see  that  within  the  second  platoon,  the  first  twelve  aircraft  need  to  be 
actuated  within  the  Oakland  Center,  whereas  the  last  eight  need  to  be  actuated  upstream 
(Salt  Lake  Center).  Since  in  general,  no  knowledge  of  the  required  boundary  conditions  is 
propagated  upstream,  we  can  predict  that  in  the  real  system,  the  last  eight  aircraft  would 
not  be  actuated  until  they  enter  the  Oakland  Center  and  that  no  solution  to  this  metering 
problem  would  be  found  without  putting  the  aircraft  on  hold. 

We  verify  this  by  simulating  this  flow.  In  Figure  14,  we  see  that  for  the  last  eight 
aircraft,  the  first  activation  time  in  the  simulator  is  higher  than  the  predicted  (upper  plot): 
the  simulated  controller  is  not  able  to  actuate  the  aircraft  on  time,  because  they  are  not  in 
its  airspace.  We  illustrate  in  the  middle  plot  that  these  aircraft  are  breaching  the  boundary 
conditions  (by  about  one  minute  each),  which  can  also  be  seen  in  the  lower  plot:  their 
delays  (the  inverse  of  the  cumulated  breaches)  become  negative,  i.e.  they  arrive  in  advance. 
This  is  an  illustration  of  distributed  and  decentralized  control:  the  actuation  occurs  in 
different  sectors,  and  the  only  communication  between  the  sectors  is  through  the  metering 

4  Small-scale  maneuvers  are  less  likely  to  be  executed  correctly  by  the  simulator  because  the  probability 
of  selecting  the  respective  aircraft  at  exactly  the  ‘right’  time  is  small,  which  explains  the  small  discrepancy 
between  the  results.  Also,  even  if  the  maneuver  is  executed  correctly  by  the  simulator,  the  resulting  flight 
plan  will  look  different  from  the  ETMS  data,  since  the  simulator  is  restricted  to  a  single  angle  of  deviation 
(0  =  22.5°). 

5The  terminology  “n  miles-in-trail”  is  a  standard  term  used  by  ATC.  It  means  that  aircraft  follow  each 
other  separated  by  n  miles. 
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conditions.  Obviously,  the  lack  of  centralized  actuation  (here  communication  and  strategic 
TMU  planning)  disables  efficient  flow  scheduling. 


3.4  Decorrelation  of  flow  metering  and  conflict  resolution 

A  working  assumption  which  is  often  made  in  flow  scheduling  is  the  decorrelation  of  meter¬ 
ing  conditions  and  cross  traffic  when  traffic  density  is  relatively  low,  that  is,  that  conflict 
avoidance  actuation  does  not  impact  metering  actuation.  We  demonstrate  this  property  of 
the  flow  by  simulating  two  streams  of  intersecting  (and  conflicting)  aircraft  (see  Figure  16). 
Two  platoons  of  aircraft  intersect  at  a  navaid  (Clovis,  in  Sector  15),  and  we  consider  pairwise 
conflicts.  One  platoon  is  subject  to  outflow  metering  conditions  (at  the  boundary  of  sector 
15  and  34),  and  we  would  like  to  quantify  the  impact  of  the  other  platoon  on  the  travel  time 
of  the  first  platoon  through  the  sector. 

We  investigate  16  different  configurations  by  selecting  the  same  inflow  conditions  for  the 
two  platoons:  ALm  G  {15,20,25,30}  milcs-in-trail,  and  the  following  outflow  conditions  for 
the  platoon  heading  towards  sector  34,  ALout  e  {15,  20, 25,  30}  milcs-in-trail  (see  Figure  16). 
We  compare  the  travel  time  for  each  platoon  without  the  presence  of  the  other,  with  the 
travel  time  when  the  two  platoons  intersect.  For  each  configuration,  we  perform  10  runs 
where  the  initial  position  of  the  aircraft  within  each  platoon  is  perturbed  by  a  uniform  noise 
of  amplitude  2  nm.  This  scenario  represents  ths  situation  in  which  each  aircraft  from  one 
platoon  conflicts  pairwise  with  an  aircraft  from  the  second  platoon.  This  results  in  160  runs.6 

For  each  (ALin,  ALout),  we  compute  the  three  following  quantities 


rji  no  cross  flow 
u  travel  time 

T'T’cross  flow 
u  1  travel  time 


=  ^Eirfe 

=  k 


■  aircraft  i 

pcross  flow,  B.C. 
'  aircraft  i 


pno  cross  flow,  no  B.C. 
■  aircraft  i 


p no  cross  flow,  no  B.C. 
'  aircraft  i 


-  due  to  cross  flow 


cross  flow  _  ro“>no  cross  flow 

u  1  travel  time  travel  time 


Here,  N  is  the  total  number  of  aircraft  ( N  =  20,  we  have  two  platoons  of  10  aircraft). 
In  each  formula  above,  ow’  B'C'  represents  the  travel  time  of  aircraft  i  in  absence 

of  the  other  platoon,  and  with  outflow  boundary  conditions  as  shown  in  Figure  16  (similar 
definitions  apply  for  the  other  variables  in  formula  (11)  and  are  obvious  from  notation). 
The  results  for  the  mean  difference  in  travel  time  are  shown  in  Figure  17  (averaged  over  10 
runs  for  each  case).  We  show  the  numerical  results  obtained  respectively  in  the  ALin  =  15 
miles-in-trail  inflow  case: 


6The  settings  for  these  runs  are:  M  £  {0.8,  0.85,  0.89},  the  vector-for-spacing  maneuver  was  limited  to 
a  5  nm  deviation  from  the  original  flight  plan.  These  settings  were  chosen  to  guarantee  short  flight  times. 
The  interval  between  controller  activation  was  set  to  A Tact  =  20sec. 
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and  for  the  AL-m  =  20  miles-in-trail  inflow  case: 


ALout 

jj'T'no  cross  flow 
u  1  travel  time 

jr^ycross  How 
u  1  travel  time 

AZdue  to 

cross  flow 

15nm 

Osec. 

28.1sec. 

28.1sec. 

20nm 

0.7sec. 

35.9sec. 

35.2sec. 

25nm 

74.2sec. 

101.9sec. 

27.6sec. 

30nm 

99.5sec. 

123.9sec. 

24.3sec. 

ATdueto  is  represented  in  Figure  17  for  the  complete  set  of  (ALin,  ALout)  investigated  here. 
Even  if  the  peak  of  ATAe  to  happens  for  (ALin,  ALout)  =  (15, 15),  the  maximum 
happens  as  expected  for  (ALin,  ALout)  =  (15,30),  which  is  the  maximal  inflow,  minimal 
outflow,  as  can  be  seen  in  (12).  The  comparison  between  (12)  and  (13)  and  Figure  17  clearly 
shows  the  predominance  of  conflict  resolution  over  boundary  conditions  for  high  density  of 
traffic  (see  last  column  in  (12)  and  (13)). 

We  see  that  the  difference  in  delay  between  separate  and  simultaneous  flow  is  significantly 
larger  if  the  aircraft  are  spaced  at  15  nm  when  compared  to  platoons  with  larger  spacing 
(see  Figure  17).  While  the  mean  difference  per  plane  is  always  larger  than  60  seconds  for 
the  15  nm-platoons,  it  is  always  smaller  than  60  seconds  for  platoons  with  larger  spacing. 
With  an  average  flight  time  of  660  seconds  over  all  scenarios,  60  seconds  corresponds  to  an 
average  delay  of  9%  in  flight  time.  The  worst  case  difference  (15  nm  inflow,  30  nm  outflow) 
is  more  than  21%  of  the  overall  flight  time.  These  numbers  are  significant,  especially  when 
considering  the  possibility  of  multiple  intersecting  platoons. 

Note  that  in  the  left  plot  of  Figure  17,  we  would  intuitively  expect  the  largest  difference 
in  delay  (cross  flow  vs.  no  cross  flow)  to  happen  for  (ALin,ALout)  =  (15,30),  which  is 
the  hardest  situation  to  achieve  (maximum  inflow,  most  restricted  outflow).  In  fact,  this 
maximum  occurs  at  (ALin,  ALout)  =  (15, 15),  which  does  not  seem  a  priori  to  be  problematic 
(conservation  of  flow).  This  can  be  explained  by  looking  at  the  middle  and  right  plot  of  Figure 
17.  In  the  absence  of  cross  flow,  the  delay  accumulated  due  to  the  boundary  conditions  is 
maximal  for  (ALin,  ALout)  =  (15,  30),  as  expected  (right  plot  of  Figure  17). When  subtracted 
from  the  difference  in  travel  time  with  cross  flow  (middle  plot  in  Figure  17),  a  max  appears 
at  (ALin,ALout)  =  (15,15),  because  5^%XtTmfle°w|(ALin,ALout)=(i5.i5)  =  0.  Thus,  these  results 
are  consistent.  The  flows  of  Figure  16  are  worst  cases,  and  it  is  rare  to  observe  such  a  density 
in  both  directions.  However,  sectors  such  as  33  are  often  subject  to  less  structured,  but  just 
as  dense  cross  flows,  and  our  predictions  are  relevant  for  this  type  of  airspace. 

Currently,  TMU  does  not  take  the  influence  of  conflict  resolution  into  account  when  mak¬ 
ing  decisions,  since  the  impact  of  ATC  on  the  sector  level  is  considered  to  be  negligible  for 
overall  flight  times.  In  the  current  system  this  is  not  always  true,  thus  leading  to  inaccu¬ 
racies  in  the  predictions  of  sector  occupancy,  i.e.  how  many  aircraft  will  be  in  a  sector  in 
t  +  15  minutes.  We  have  thus  shown  in  a  particular  case  that  the  influence  of  conflict  resolu¬ 
tion  actuation  increases  with  higher  traffic  density,  which  therefore  necessitates  information 
feedback  from  the  sector  to  the  Center. 
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4  Conclusion  and  Current  Work 


We  have  derived  a  control  theoretic  model  of  sector-based  traffic  flow  using  hybrid  automata 
theory.  We  used  a  subset  of  this  model  to  generate  Lagrangian  analytic  predictions  of  the 
traffic  flow  (dynamic  sector  capacity,  extend  of  traffic  jams),  which  were  linked  to  Eulerian 
models  of  the  NAS.  We  used  these  results  to  predict  the  conditions  under  which  airspace 
saturation  cannot  be  treated  at  the  level  of  a  single  sector,  but  requires  centralized  actuation 
(communication  and  strategic  TMU  planning).  These  predictions  were  validated  against  an 
abstraction  of  the  real  system  (a  simulator  using  the  full  model,  which  we  validated  against 
real  data).  Finally,  we  generated  flow  conditions  under  which  we  can  decorrelate  metering 
from  conflict  resolution. 

This  last  result  is  useful  for  our  current  work,  because  it  enables  us  to  ignore  conflict 
resolution  for  scheduling.  Our  current  work  [23]  aims  at  generating  fast  algorithms  to  design 
metering  protocols  for  this  converging  airspace.  The  emphasis  of  this  work  is  in  design, 
rather  than  modeling,  and  it  uses  all  modes  of  the  automaton  shown  in  Figure  4.  Therefore, 
no  analytical  solution  is  available  to  solve  this  problem.  In  fact,  the  linear  program  derived 
in  section  2  generalizes  to  a  Mixed  Integer  Linear  Program  (MILP),  which  we  believe  to  be 
NP-complete  in  the  general  case.  However,  in  specific  cases  relevant  to  ATC  automation,  we 
have  been  able  to  derive  polynomial  time  algorithms  to  solve  this  program,  which  leads  us  to 
believe  that  this  method  could  be  applied  to  real  time  scheduling  and  airspace  decongestion. 
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Figure  1:  Control  hierarchy  in  the  current  structure  of  NAS. 
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350 


Figure  2:  ATC  sectors  modeled  for  this  study:  32,  33,  34,  13  and  15  within  the  Oakland  ARTCC  (labeled 
above  as  ZOA32,  etc.).  Most  crucial  jetways  and  waypoints  for  a  San  Francisco  approach  are  shown  here. 
The  data  modeled  comes  from  FACET  [11]  as  well  as  JEPPESEN  high  altitude  en  route  charts  [24]. 
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Figure  3:  Visual  display  of  the  simulator.  Traffic  in  the  Oakland  ARTCC  (aircraft  not  in  this  Center  have 
been  filtered  out).  This  plot  has  been  generated  using  ETMS  data. 
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Figure  4:  Hybrid  automaton  representing  the  action  of  one  controller  on  a  single  aircraft.  Each  of  the  eight 
modes  represents  one  possible  state  of  the  aircraft.  The  arrows  joining  these  states  are  the  mode  switches, 
initiated  by  the  controller.  The  dash-dotted  transitions  are  used  for  the  analytical  solution.  The  complete 
set  of  arrows  is  used  for  the  simulation. 
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Figure  5:  Overlay  of  trajectories  merging  into  San  Francisco  (11  hours  of  traffic).  The  data  modeled  comes 
from  ETMS  and  FACET  [11].  Our  Lagrangian  approach  models  each  of  the  trajectories  included  in  this 
plot. 
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Figure  6:  ATC  control  action  on  the  merging  flow.  The  traffic  jam  extends  from  SFO  to  the  edge  of  the 
platoon  with  ATout  •  iWn  miles- in-trail  (one  aircraft  every  A Tout  time  units) .  Once  the  actuation  point  (edge 
of  the  traffic  jam)  has  moved  upstream  (into  sector  33),  sector  34  is  called  saturated. 


Figure  7:  Shock  construction.  The  aircraft  trajectories  are  represented  in  the  (x,t)  plane.  They  originate 
at  t  =  0  from  the  horizontal  axis  (white  circle  on  each  trajectory).  After  some  amount  of  time,  the  aircraft 
may  be  switched  to  speed  vmin  at  location  (x®wltch,  f®wltch)  (shaded  circle  on  each  trajectory).  Ultimately 
they  reach  xex>  the  entrance  of  TRACON  (black  circle). 
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Figure  8:  Switching  curve  (shock)  for  a  vanishing  traffic  jam.  x  denotes  the  distance  to  the  metering 
point  (SFO).  The  lines  are  the  trajectories  of  the  aircraft  in  the  (x,t)  space.  The  positions  of  aircraft  are 
represented  every  1000  sec.  as  dots.  Once  they  have  passed  through  the  shock,  they  are  separated  by 
i!mjr,  ATout .  The  point  {xm,tm)  is  the  furthest  reachable  point  by  this  traffic  jam.  Note  that  the  slope  of  the 
lines  changes  through  the  shock.  The  slope  difference  can  hardly  be  seen  visually,  because  the  speed  change 
is  small. 
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Figure  9:  Program  flow  of  the  simulator. 
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Aircraft-Cost  Landscape 


Figure  10:  Top:  cost  values  for  all  possible  maneuver  combinations  in  a  two-aircraft  intersection  scenario, 
where  the  eight  maneuvers  of  Figure  4  are  enabled  (thus  generating  82  =  64  possible  values  of  J).  Four  out 
of  64  examples  are  extracted  and  illustrated  (four  lower  pictures),  (a)  Both  aircraft  maintain  same  speed; 
(b)  Aircraft  A  takes  a  shortcut  maintaining  aircraft  B  at  max  speed;  (c)  A  makes  a  VFS  at  low  speed;  (d) 
A  does  nothing,  B  is  not  able  to  prevent  the  loss  of  separation.  In  this  case,  the  simulated  controller  would 
choose  solution  (b)  since  the  lowest  cost  is  associated  with  that  maneuver.  The  cost  J  has  been  truncated 
at  5  •  103  for  readability.  Thus  a  LOS  cannot  be  visually  differentiated  from  a  breach  in  this  plot,  though  it 
can  in  our  data.  Thus,  too  large  breaches  as  well  as  LOS  do  not  appear  on  this  plot. 


Figure  11:  Flight  time  comparisons  for  the  first  100  aircraft  going  through  sector  33  in  the  ETMS  data  set 
we  used.  The  dots  are  the  flight  times  for  the  ETMS  recorded  points.  The  solid  curve  is  the  result  of  our 
simulations. 
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Figure  12:  An  example  of  maneuver  caused  by  conflict  resolution,  reproduced  by  the  simulator.  The 
recorded  data  (dashed)  exhibits  an  actual  shortcut  from  the  filed  flight  plan  (solid).  The  simulated  trajectory 
(dashed-dotted)  is  a  shortcut  of  the  same  type. 
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Figure  13:  Sector  33:  traffic  flow  for  the  merging  traffic  simulation  of  Figures  15  and  14.  The  radius  around 
the  aircraft  is  2.5  nm.  The  solid  lines  represent  the  aircrafts’  flight  plan.  The  dotted  lines  correspond  to 
maneuvers  assigned  by  the  simulator.  The  simulator  makes  extensive  use  of  the  vector-for-spacing  to  meter 
the  aircraft. 
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Figure  14:  First  actuation  times  of  the  aircraft  (upper  plot),  simulated  and  predicted;  breaches  in  metering 
conditions  (middle  plot),  simulated;  delays  (lower  plot),  simulated,  for  the  case  of  the  two  platoons  of 
Figure  15. 
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Figure  15:  Shocks  generated  by  two  successive  platoons.  The  first  shock  is  steady  in  time  (it  only  propagates 
backward  in  space).  It  corresponds  to  a  piling  up  process  on  a  highway  where  all  vehicles  slow  down  at  the 
same  time.  The  second  shock  propagates  backward  in  space  and  time  (which  is  much  harder  to  handle  in 
practice,  because  actuation  must  be  performed  upstream  first). 
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Figure  16:  Two  intersecting  platoons  at  the  Clovis  navaid  for  the  configuration  in  which  the  aircraft  have 
an  inflow  of  20  miles-in-trail.  Aircraft  are  pairwise  conflicting  (aircraft  i  from  platoon  1  conflicts  with  aircraft 
i  of  platoon  2),  assuming  the  flight  plan  is  to  follow  the  straight  line  at  the  initial  speed.  The  simulator  is 
thus  forced  to  adjust  each  of  the  pairs.  Deviations  (vector-for-spacing)  are  barely  visible  on  the  plot  because 
of  their  small  amplitude  (5nm). 
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Figure  17:  Left:  Plot  of  ATdue  to  (in  sec.)  averaged  over  10  runs.  Middle:  Plot  of  <5Ttcr^f,  (in  sec.) 

cross  flow 

averaged  over  10  runs.  Right:  Plot  of  ^^TraveUime™  (in  sec-)  averaged  over  10  runs. 
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